Entanglement of Solitons in the Frenkel-Kontorova Model 
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We investigate entanglement of solitons in the continuum-limit of the nonlinear Frenkel-Kontorova 
chain. We find that the entanglement of solitons manifests particle-like behavior as they are char- 
acterized by localization of entanglement. The von-Neumann entropy of solitons mixes critical 
with noncritical behaviors. Inside the core of the soliton the logarithmic increase of the entropy is 
faster than the universal increase of a critical field, whereas outside the core the entropy decreases 
and saturates the constant value of the corresponding massive noncritical field. In addition, two 
solitons manifest long-range entanglement that decreases with the separation of the solitons more 
slowly than the universal decrease of the critical field. Interestingly, in the noncritical regime of 
the Frenkel-Kontorova model, entanglement can even increase with the separation of the solitons. 
We show that most of the entanglement of the so-called internal modes of the solitons is saturated 
by local degrees of freedom inside the core, and therefore we suggest using the internal modes as 
carriers of quantum information. 



I. INTRODUCTION 

It is well known that nonlinear theories manifest a va- 
riety of classical configurations such as solitons and in- 
stantons and dynamical effects such as symmetry break- 
ing, confinement and phase transitions, which can never 
be seen using an ordinary perturbative expansion. Such 
exact classical configurations can provide us with pow- 
erful insight into the nonperturbative content of those 
theories and may exhibit novel emergent particle-like be- 
havior. Classical configurations such as solitons are com- 
monly observed solutions of nonlinear relativistic or non- 
relativistic equations. They are localized, finite energy 
solutions, which can travel while retaining their shape. 
These classical configurations also strongly determine 
the structure of the corresponding quantized theories. 
For example, the Hilbert space of the well-known one- 
dimensional 4 and sine-Gordon field theory models are 
classified to topologically inequivalent sectors, which can- 
not be obtained by power expansion of the global vacuum 
state. However, these sectors can be obtained by quan- 
tizing a field around the nontrivial backgrounds of either 
the <p 4 kinks or sine- Gordon solitons 

It is interesting to ask whether nonperturbative effects, 
such as collective particle- like behavior, arc manifested in 
the entanglement structure of nonlinear theories. Much 
of the recent study of entanglement in multipartite sys- 
tems has been carried out in linear systems such as har- 
monic and spin chains @, 0j- In such one-dimensional 
system it has been found that entanglement is a measure 
of criticality and manifests a universal scaling behavior 
[3, [Ej- The behavior of entanglement in nonlinear sys- 
tems may challenge some of the insights from the linear 
systems, and shed new light regarding the nature of en- 
tanglement. The study of entanglement in nonlinear sys- 
tems may also suggest new implementations of quantum 
information processes. 

In this work we investigate entanglement in a chain of 
nearest neighbors interacting oscillators subjected to a 
periodic on-site external potential, known as the Frenkel- 



Kontorova (FK) model d, 0]- This model possesses rich 
behavior and so we parameterize it from two aspects: 
its weak/strong-coupling regime and critical/noncritical 
regime. The continuum-limit solutions of the FK model 
are the well-known sine-Gordon solitons. Although we 
deal with a finite discrete system, it still retains most of 
the characteristics of the continuum and therefore in the 
following we shall use the same field-theory terminology 
(vacuum-sector, soliton-sector, etc.). 

Roughly speaking, one can say that a soliton solution 
is one in which the effective (mass) 2 of the quantum free- 
field around the classical background decreases inside the 
core of the soliton and even becomes negative. Outside 
the soliton the field behaves as a regular massive free 
theory. Therefore, the field interpolates critical-like be- 
havior inside the soliton's core and noncritical behavior 
outside the core. 

It can be expected, therefore, that the von-Neumann 
entropy of a region which coincides with the soliton will 
manifest a mixture of critical and noncritical behaviors, 
depending on its size. This is the first result we pro- 
vide. We compute the entanglement in the single-soliton 
sector ground state. We find that as long as the block 
is confined to the soliton's core, entanglement increases 
logarithmically but faster than the critical harmonic sys- 
tems. However, as we keep increasing the size of the 
block so that it grows outside the soliton's kink, the en- 
tanglement decreases and saturates the asymptotic con- 
stant limit that corresponds to the noncritical massive 
system. Therefore, the soliton manifests a particle-like 
localization of entanglement. 

Next we study long-range effects of entanglement, that 
is the entanglement between two spatially separated 
blocks, and use the logarithmic negativity as a measure 
of entanglement. We find that the entanglement distri- 
bution is peaked around the location of the solitons. We 
then study soliton-soliton entanglement as a function of 
their separation and compare the results with the cor- 
responding behavior in the vacuum sector. The latter 
is also equivalent to the ground state of a linear har- 
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monic chain with mass m ~ 1^/5, where g is the coupling 
strength. 

In the two-soliton sector we see complex behavior of 
long-range entanglement. We find that the entanglement 
between two solitons is larger than the corresponding en- 
tanglement in the vacuum sector. For relatively small 
separations of the solitons the entanglement decays with 
the separation more slowly than that of the critical vac- 
uum sector. Moreover, in the noncritical limit entan- 
glement loses one of its most fundamental properties: it 
is no longer a monotonically decreasing function of the 
distance between the solitons! 

Much of our discussion is devoted to the strong- 
coupling limit. In the weak-coupling limit, where the 
kinks become impurities that are confined to a few par- 
ticles, we find that the entanglement loses both its local- 
ization and long-range effects. 

Finally, we suggest and investigate the possible use of 
solitons as carriers of quantum information. We propose 
utilizing the internal localized modes, which describe col- 
lective and vibrational perturbations of the classical soli- 
ton, as the carriers of quantum information. The soliton 
is handled as a non-perturbative classical object while 
the internal modes are quantized linear perturbations 
"around" it. Since the internal modes are "attached" 
to the soliton, by moving the soliton, quantum informa- 
tion can be transported in the system. (Transportation 
of quantum information has been discussed for harmonic 
chains fl] , spin chains and Josephson arrays [lO] ) . Wo 
also formulate multi-soliton processes, which generate en- 
tanglement. 

The paper is organized as follows. In the next section 
we describe the FK model and analyze its strong/weak- 
coupling and critical/noncritical limits. Then in section 
3 we investigate localization of entanglement in the con- 
tinuum limit of the FK model. In section 4 we present 
the soliton-soliton long-range entanglement. In section 5 
we briefly discuss the loss of both long-range effects and 
localization in the weak-coupling limit. In section 6 we 
discuss possible quantum information implementations. 
Lastly, in section 7 we summarize the significant results 
and outline future directions of research. 



II. FRENKEL-KONTOROVA MODEL 

We consider N particles with canonical coordinates <f> n 
and 7r„, described by the Hamiltonian of the FK model 
[1,0: 

N 1 

H = E O^n + 2V sub (<f> n ) +g{4> n+ l - Cbn) 2 ), (1) 



where <p n is the displacement of the nth particle from its 
equilibrium position, g is the coupling constant between 
particles in the chain and 



We use dimensionless units in which the particles mass 
equals 1 and the period and amplitude of the substrate 
potential V su b(<fi) are a s = 2tt and e s = 2 respectively. 
In the ground state of the system all particles occupy 
minima of the potential such that 4> n = or <f> n = 2ir, etc. 
Then the simplest excited state is a soliton that connects 
two neighboring ground state solutions as illustrated in 
figure [TJ 



A. Classical solutions 

The exact solutions of the discrete FK chain are not 
known in an explicit analytical form due to the discrete- 
ness of the problem. However, in the strong coupling 
limit, g » 1, an approximate discrete solution can 
be obtained using the continuum limit approximation, 
n — > x — na s and <f>i(t) — > (j>(x,t) [HI]. In this framework 
we begin with the sine-Gordon Hamiltonian density 

H(x,t) = ln 2 +ga 2 s ±(<f>>) 2 + V sub (cb), (3) 

where a dot or a prime represents differentiation with 
respect to time or space respectively. The equation of 
motion is then the integrable sine-Gordon (SG) equation: 



d 2 (j) 



ga s 



dx 2 



sin d> = 0. 



(4) 
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FIG. 1: (color online). Schematic structure of the FK chain 
for TV = 10. The blue circles are the vacuum equilibrium loca- 
tions, where the black squares represent the soliton solution. 



Imposing the boundary conditions 



dx 



(±oo) = 



(5) 



and rescaling x — > x/a s , we obtain the SG soliton solu- 
tion: 



/ x 

<t>SG{%) = 4 tan -1 exp I — a— 



X 



V9 



(6) 



V su b(<j>) = 1 - cos(/>. 



(2) 



where ^fg defines the soliton's width in units of the pe- 
riod of the substrate potential, a = ±1 and X is the 
coordinate of the soliton's center. For g > g con ~ 16 the 
continuum limit is a reasonable approximation. We can 
then describe our iV-particle system using the continuum 
solution by sampling ipsa at N points with separation 1, 
where X is sampled at the middle. Throughout the paper 
we shall refer to this sampled classical solution (Eq. [5]) 
simply as the soliton solution (for g > g CO n)- Note that as 
we increase g, even though the separation between par- 
ticles is kept constant, the number of particles (points) 
that sample the core of the soliton, y/g, increases. 
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FIG. 2: (color online). Single soliton (fisG and double-soliton 
4>2SG solutions of the finite SG equation sampled at N = 1000 
points. The length of the region is given by 2L = 16Aj where 
Xj is the Josephson length. For the single soliton solution the 
stability region for the magnetic field is H 6 [0.0019, 1.0052) 
and for the double-soliton H € [0.1023, 1.0622). See details 
in Appendix A. 



In order to investigate long-range entanglement be- 
tween solitons, we are interested in exploring double- 
soliton solutions as well. There are no multi-soliton static 
solutions in infinite systems. However, the finite sine- 
Gordon equation also possesses multi-soliton static solu- 
tions (Eq. [4]). This equation describes the magnetic flux 
in a long Josephson junction of length 2L with a constant 
homogeneous external magnetic field H, perpendicular to 
the barrier and with bias current line density I. We use 
these double-soliton solutions, 4>2SG, as the continuum 
limit approximation of the of FK model in a finite sys- 
tem. In appendix A we describe the analytical solutions 
of the double-soliton solutions. For illustration purposes, 
we present in figure [2] single soliton and double-soliton 
configurations. More details are provided in Appendix 
A. 

In the weak-coupling limit of the FK model (g < g con ) 
there are no known analytical solutions. Still, in this 
regime there are kink solutions, 4>kink, which are local- 
ized to a small number of particles. These kink solutions 
can be found numerically. We use a descending gradient 
algorithm that minimizes the classical configuration en- 
ergy, given an initial solution to start with. Apparently, 
in the weak-coupling limit, where g is small (g ~ 1) the 
analytical continuum term is still a good enough starter 
for the minimizer to find the minimal stable configura- 
tion. 



B. Quantization 

In order to explore quantum mechanical effects in soli- 
tons we use a semi-classical framework. For simplicity, 
we provide the analysis on the continuous solutions. We 
consider small perturbations rj(x, t) on the background of 



the classical solution <j>o (x) (single soliton, double-soliton 
or kink): <f>[x) = <j) a {x) + r](x,t) Q. Setting H = U + T, 
where U is the potential energy and T is the kinetic en- 
ergy, we assume the solutions are static, T = 0. Then we 
can make a functional Taylor expansion of the potential 
energy 



U(tf>) 



dx 



WY + VsubW 



about 4>q: 



dx- 



{ v{x ){- g ^ + (^)^{x)}, (7) 



where cubic and higher terms in r\ are neglected. This 
requires that the magnitude of r\(x) as well as the third 
and higher derivatives of V(<p) at (f>Q be small. 

The eigenvalues and eigenfunctions of the operator 
(— gV 2 + d 2 V su b / d(j) 2 ) evaluated at <j>(x) = <j>a(x) are then 
the generalized solutions of the Schrodinger-like equation 



SV 2 



d 2 V s , 



<po(x) 



rji(x) =ufrji(x), (8) 



where the rji (x) are the orthonormal normal-modes of the 
fluctuations around <fio(x). In the continuous solution I 
is a continuous index with possible additional discrete 
values. 

Translating the above analysis to the discrete system 
we obtain a set of equations, that is, the rji(n) and ujf 
are obtained by diagonalizing the matrix B(n, n) = 2g + 



cos</>o(n), B(n, nil) 



-g. We require that wf > 0; 



otherwise the system is in a nonstable configuration. 




FIG. 3: (color online). Normal modes in the background of 
the single soliton solution, shown in figure [2] The four first 
normal modes and the 10'th normal mode are shown. The 
first normal mode is an internal mode. 

We obtain the quantum-mechanical modes of the sys- 
tem by setting: 



(9) 
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FIG. 4: (color online). Normal modes in the background of 
the double-soliton solution shown in figure [2] The four first 
normal modes and the tenth normal mode are shown. The 
first two normal modes are internal modes. 



and quantizing the normal- mode coefficients c/. The 
canonical degrees of freedom are then: 



1) = Mn) + ^^H T7 ^=(a i e^ t + 



7r(n) 




a, e 1 



), 



(10) 



a\e~ luJlt 



where a\ and a\ are annihilation and creation operators 
of the normal modes, and 



E 



U{<h) 



1 



(ii) 



where n/, is the excitation number of the Z'th normal 
mode. 

The background classical solution (j>o could either be 
an absolute minimum or a local minimum of U ((f)). Con- 
sider first the case that 0o is an absolute minimum. Then 
all particles are located at their equilibrium positions and 
the solution does not contain kinks or solitons. The nor- 
mals modes are then identical with the "phonon modes" 
of a harmonic chain: 

rji(n) (x exp(jfc;n); ki — 2ttI/N (12) 

and the spectrum is given by [l2j : 

ojf = l + 2s-(l -cos h). (13) 

Upon quantization, the quantum states of the FK sys- 
tem in the absolute minimum case, which we also refer to 
as the " vacuum sector" , coincide with that of a harmonic 
chain. Hence, in the continuum limit the vacuum sector 
of the FK model corresponds to a free massive scalar 
field with m ~ 1/ ^J~g. and the ground state has the same 
properties of the corresponding free field vacuum state. 

If 4>o is a local minimum, then the solution contains 
solitons (or kinks). In the continuum limit the eigenvalue 
equation is 



92 

9 Q^2 + C °^SG 



r)i(x) = ujfm(x)- 



(14) 



In the limit of an infinite system the translational invari- 
ance of the system gives rise to a zero frequency mode 
(Goldstone mode, or "zero mode"), which requires spe- 
cial attention. This problem is avoided in our case due to 
discreteness and finiteness of the system. Nevertheless, 
in correspondence with the infinite limit, the zero mode 
is still characterized by a localized shape of a bound state 
and will be referred to as an internal mode. (In other non- 
linear models which have soliton solutions such as the (j) 4 
model, there can be additional internal localized modes 
which are not zero modes.) 

For clarification, we present several normal modes and 
their eigenfrequencies for the solutions presented in figure 
[2] In figure[3]the first 4 normal modes and the 10'th nor- 
mal mode of the single-soliton solution are shown, where 
the first 4 eigenfrequencies are lo\ ~ 0.0007, 0J2 — 1.0002, 
~ 1.05, UJ4 ~ 1.2, given in dimensionless units. Note 
that the first eigenvalue is less than 1 and corresponds to 
the internal mode, while all the others are phonon modes. 
In figure 2J the first 4 normal modes and the 10'th nor- 
mal mode of the double-soliton solution are shown, where 
the first 5 eigenfrequencies are u>\ ~ 0.008, 102 — 0.016, 
ui3 ~ 1.013, LU4 ~ 1.07, ~ 1.3. Here the first two 
normal modes are internal modes. 

As we shall see, it is the presence of the internal modes 
which gives rise to qualitatively different behavior of the 
quantum mechanical states in the solitonic sector, such as 
particle- like behavior and localization of entanglement. 



C. Correlations in the vacuum and soliton sectors 

We would like to compare the behavior of correlations 
in the ground states of the vacuum and soliton sectors in 
the strong and weak coupling regimes. Let us define the 
correlation lengths £ n = (4>o4> n ) and v n — (^0^71)1 where 
n is the (integer) distance between the particles. 

We begin with the strong-coupling regime. In the vac- 
uum sector the correlation lengths , v n have been thor- 
oughly analyzed in connection with the harmonic chain 
These correlations are classified by the critical and 
noncritical regimes. In the critical limit there are long- 
range correlations as there is no length scale in the system 
(to — > 0). In this limit the correlations scale as: 



'log(-), 
n 

1 



(15) 



In a system with N particles and a coupling strength g 
the critical regime is valid if N << ^fg. 

The noncritical regime is characterized by a finite 
length scale. Here the correlations decay exponentially 
with the distance: 



A/2 ' 



(16) 



,3/2' 
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The noncritical regime is defined for N » Jg. For 
N = 1000 [13], for example, the system is critical for 
g > 10 7 , and noncritical for g con < g < 10 5 (l4j |. 

We turn now to the soliton sector. For the critical 
limit the number of particles in the core of the soliton 
should be greater than the number of particles in the 
whole chain: N < ^/g. Therefore, in the critical limit 
the vacuum sector and the soliton sector coincide: for 
a constant N taking g to infinity causes the loss of the 
internal mode, which now enters into the phonon band. 
We define g m ax(N) as the maximal coupling constant for 
which there is a discrete internal mode. For g > g max (N) 
no difference is predicted between the ground states of 
the vacuum and soliton sectors since their spectra and 
normal modes become practically identical. 

On the other hand, in the weak-coupling limit, g < 
g C on, the correlations decay very fast and again one may 
expect no significant difference between the vacuum sec- 
tor and the kink sector. Therefore, we expect qualitative 
difference between the vacuum sector and the soliton sec- 
tor only in the noncritical, strong coupling regime. 
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FIG. 5: (color online). Log-log scale of f n , where < n < 500. 
The plots are given for the vacuum sector and the soliton sec- 
tor in the critical (g = 10 10 ), noncritical (g — 10 4 ) and weak- 
coupling (g = 5) regimes. In the soliton sector the center of 
the soliton is located at n — 0. 

In figure we show on log-log scale £„ for N — 1000, 
where < n < 500. The plots are given for the vacuum 
sector and the soliton/kink sector of the FK model, where 
the center of the soliton/kink is located at n = 0. We 
study the following coupling constants: g — 10 10 (critical 
limit), g = 10 4 (noncritical regime) and g = 5 (weak- 
coupling limit). First of all, we observe that the vac- 
uum and soliton/kink sectors are indeed undistinguish- 
able both in the critical limit and the weak-coupling 
limit. In the noncritical regime, however, £* oi is qual- 
itatively different from £™ c . Inside the soliton's core 
(n < ^/g = 100), I;™ 1 decays more slowly than £™ c and 
even ^ ltlcal ^ As n becomes greater than ^/g, however, 
we see a strong decay in , which quickly becomes par- 
allel to the graph of £™ c in the noncritical regime. 



We explain this behavior on physical grounds by re- 
lating to the effective (mass) 2 of the linearized free-field 
around the critical vacuum and the (noncritical) soliton 
sector: 



'"'critical u > 

™lol(. x ) ~ COS(0oO)) = 



1 X = 

7f. \A»y/9 



(17) 



As (mass) 2 become negative for x < yfg the correlation 
length inside the core decays more slowly than that of the 
critical vacuum sector. Outside the core the correlation 
restores the massive vacuum sector behavior. We expect 
that the entanglement of solitons will have a similar be- 
havior to the correlation length. 



III. LOCALIZATION OF SOLITON 
ENTANGLEMENT 

We turn now to the entanglement of the FK model in 
the strong-coupling limit. Both in the vacuum sector and 
in the soliton sector we assume the normal modes are in 
their ground states, therefore the quantum-mechanical 
states are Gaussian. Following let us represent the 
local canonical variables of our N-mode system (Eq. [T0|) 
by the vector 

y = (<M) T , 

where <\> = (<fo, <f> 2 , . . . , 4> N ) and n = (m,^, ■ ■ ■ ,n N ). 
The commutation relations may thus be expressed as 

[Va,Vp] =iJ a p, 
where J is the so-called symplectic matrix: 



J = 



1 
-1 



Assuming (y) — 0, the state of the system is entirely 
characterized by the matrix of the second moments, the 
so-called phase-space 2N x 2N covariance matrix (CM): 

M(y)=Re( yy T ). 

Under a symplectic transformation y — Sy a Gaussian 
state M is mapped into a Gaussian state M = SMS T , 
where SJS T — J. A theorem due to Williamson [TBI, 
[l6j states that there always exists a certain symplectic 
transformation Sw that brings M to the normal form 
("Williamson form"): 



W = S W MS^ = diag(Ai, A 2 



, An, Ai, A2, . . . , An), 



where the diagonal elements Ai, A2, . . . , Ajv are referred as 
the symplectic eigenvalues and must be greater or equal 
to 1/2 according to the uncertainty principle. 

Now, suppose our pure iV-mode system (Eq. [10]) is 
partitioned into two sets y& and ys- In order to mea- 
sure the entanglement between parts A and B we bring 
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the reduced covariance matrices Ma(ua) and Mb(Ub) 
into their Williamson normal forms. Both Ma{va) an d 
Mb(ub) have the same Xj > 1/2 symplectic eigenvalues. 
The entanglement is then measured by the von-Neumann 
entropy Es'- 

E S = J2 S ^)> ( 18 ) 

\ 3 



where 




Since the <fi — it correlations vanish in our model, Xj are 
given by the square roots of the eigenvalues of HaGa (or 
G a Ha) 1, where 

G=(<W> T ), H=(7nr T ). (20) 
A. Mixing critical and noncritical behaviors 
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FIG. 7: (color online). Localization of entanglement in the 
soliton sector of the FK model, logarithmic scale. The in- 
crease of entanglement inside the soliton is faster than the 
corresponding one in the vacuum sector. Both entropies co- 
incide outside the soliton. The entanglement of the soliton 
sector has a maximal value and is therefore localized. Note 
that the logarithmic increase inside the soliton is even faster 
than the corresponding one in the critical vacuum sector. 
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FIG. 6: (color online). Localization of entanglement in the 
soliton sector of the FK model. In the noncritical regime 
(g — 10 4 ) the size of the core is ^fg = 100. The critical limit 
is given for g — 10 8 . Localization of entanglement is seen in 
the noncritical regime, where maximal entropy is obtained for 
I = V9- 

In figure[S]we plot the von-Neumann entropy of a block 
of I particles in the chain where in the soliton sector the 
center of the soliton coincide with the center of the block. 
The plots are given for both the vacuum sector and the 
soliton sector (in the strong-coupling regime of the FK 
model). The plots are given for two regimes: critical 
(g = 10 8 ), in which case both sectors have the same 
behavior and noncritical (g = 10 4 ). 

First of all, let us examine the behavior in the vacuum 
sector. In the critical regime the entanglement increases 
as the size of the block increases. This effects stops, how- 
ever, due to the finiteness of the chain (due to edge effects 
this happens long before I = N/2). In the noncritical 



regime the entanglement increases for a finite region and 
then saturates a constant value. 

In figure [JJ the same graph is plotted on logarithmic 
scale for I < 100. We see that the critical vacuum sec- 
tor fits the linear curve for I < 10 15 ~ 30, where the 
noncritical vacuum sector fits the linear curve for / < 16. 

In the soliton sector we observe in figure [6] qualita- 
tively different behavior, as anticipated from the discus- 
sion on the correlation length in the previous section. 
The von-Neumann entropy of the soliton mixes both the 
critical and noncritical scalings: it scales logarithmically 
inside the soliton, that is for I < yfg = 100, and then 
it decreases to the noncritical vacuum sector asymptotic 
value. As seen in figure [71 inside the soliton the logarith- 
mic increase of the soliton sector is faster than the cor- 
responding vacuum sector. This characteristic behavior 
implies that the entropy in the soliton sector must have a 
maximal value. In contrast to the vacuum sector, the en- 
tanglement in the soliton sector is localized as it reaches a 
maximal value at a certain length scale corresponding to 
the size of the soliton. Note that the logarithmic increase 
of the entanglement inside the core of the soliton is even 
faster than that of the critical vacuum sector. This result 
corresponds to the effective negative (mass) 2 inside the 
core. 

For reference we plot in figure [8] the von-Neumann en- 
tropy of the soliton sector for g — 200, g = 1000 and 
g = 25000. We may say that entanglement localization 
" measures" the size of the soliton, as maximal entangle- 
ment is achieved for / = ^fg in all three cases. We can 
see that as the system enters deep into the noncritical 
regime (g = 1000, 200) the entanglement becomes more 
localized and reaches a higher maximum. 
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Therefore, as g decreases Eg° l falls more rapidly outside 
the core and the initial logarithmic jump of the entan- 
glement becomes sharper. This effect stops, however, 
for g < N, as then already for 1 = 1, the entanglement is 
close to its maximal value. As anticipated from the corre- 
lations behavior, for any value of g in the strong-coupling 
regime (g > g C on), the logarithmic increase is even faster 
than the corresponding one in the critical vacuum sec- 



tor: a sol {g) > a vac (g 
is obtained for g = N. 



oo 



where a^ x (g) ~ 4/9 



C. Maximal entropy 



FIG. 8: (color online). Localization of entanglement in the 
soliton sector for g — 200, 1000, 25000. Entanglement local- 
ization "measures" the size of the soliton, as maximal entan- 
glement is achieved for / = ^fg in all three cases. 



B. Logarithmic prefactor 
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FIG. 9: (color online). The prefactor a(g) of Es = a log? 
for I << y/g. Interestingly, for any value of g in the strong- 
coupling regime {g > goon), a sol (g) > a vac (g -> co) = 1/3. In 
the soliton sector a^ l ax (g) ~ 4/9, and is obtained for 



N. 



We would like now to obtain the prefactor a(g) in the 
logarithmic scaling Es = &(g) log / for small enough val- 
ues of I, I << y/g, in the vacuum sector and the soli- 
ton sector. As a first order for the linear curve we take 
the first two points in the graph: (log 2, -Eg (log 2)) and 
(log4, £?s(log4)). For g — > oo both models have approx- 
imately the well-known universal factor of a c — 1/3 in 
critical bosonic one-dimensional fields. In figured] we plot 
a(g) for g < g max (N = 1000) ~ 1.3*10 5 , where g is given 
in a log 10 scale. We observe that as g gets closer to N 
the initial logarithmic jump of the entanglement becomes 
sharper. This is explained from the following conditions: 

, sol 
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FIG. 10: (color online). Maximal entropy as a function of g. 
uji is minimal for g — 200, exactly where the von-Neumann 
entropy reaches its maximal value. 



Maximal localized entanglement is achieved in the soli- 
ton sector for lower values of the coupling constant g 
(yet not too low, so that the strong-coupling limit still 
holds). We show that the maximal value of the entan- 
glement is determined by the frequency of the internal 
model. In figure [TOl the maximal localized von Neumann 
entropy as a function of g is shown, where for reference, 
we plot on the same graph 21og(cji) + 15, so that both 
plots appear with the same scale. Clearly, lo\ is mini- 
mal exactly where the von-Neumann entropy reaches its 
maximal value (g = 200). 

In passing we note that as g gets closer to the critical 
limit, Ldi increases and eventually penetrates the phonon 
band and stops behaving as a bound state. In addition, 
uj\ increases as discreteness increases (g decreases). This 
is associated with the existence of a potential barrier in 
a discrete lattice that must be overcome to move the 
soliton along one lattice spacing, and is known as Peierls- 
Nabarro potential [l7|. Note that u>\ has an oscillatory 
behavior for g > 200. Oscillations of 0Ji(g) are known 
for models that deviate from the FK sinusoidal potential 
(21) [3 Ell- Here we observe the oscillations in the FK model 
for large values of g, where u)\ starts to increase. 
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IV. LONG-RANGE SOLITON-SOLITON 
ENTANGLEMENT 

We now explore long-range effects of entanglement in 
the FK model. We are interested in the entanglement 
between two separated blocks of particles A and B with 
the same number of particles I, separated by d parti- 
cles. Now Mab does not constitute a pure-state and 
the von-Neumann entropy cannot be used anymore. We 
therefore use a mixed-state entan glem ent measure - the 
logarithmic negativity (LN) [2CJ, |2lJ. This measure is 
based on the following observation: when the parts are 
entangled, reversing time direction in one part of the sys- 
tem: y^B = {<t l A,4>B,'n'A,—'n'B) T , breaks the symplectic 
symmetry of their covariance matrix, yielding Xj < 1/2, 
which violates the uncertainty principle. The amount by 
which Xj < 1/2 is a measure of the entanglement: 

E LN = -Y^H2X 3 ). (22) 

3 

A. Distribution of long-range entanglement 
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FIG. 11: (color online). Long range entanglement between 
two sliding blocks in the presence of two solitons with g — 
3000, / = l so i = 56 and d = d so i — 281. LN is maximized 
when the blocks coincide with the solitons' cores (the average 
of the centers is located in the center of the chain, n — 0) . 

In the double-soliton sector of the FK model we explore 
long-range entanglement between two solitons, which are 
characterized by the length scale l so i ~ yfg. We obtain 
different values of the distance between the solitons, d so i , 
by choosing appropriate values for H (and hence k) as 
described in appendix A. 

We first show that the most significant contribution 
to long-range entanglement originates from the solitons. 
Consider a certain double-soliton configuration with soli- 
tons size l so i and separation d so i. Now let us take two 
blocks with the same size and separation: / = l so i, 
d = d so i and slide them along the chain. We expect that 
as the blocks' positions become identical to the solitons' 
positions we will get the maximal value of LN. In figure 
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-500 -375 -250 -125 

FIG. 12: (color online). First and second participation func- 
tions of the double-soliton and vacuum sectors, shown in log- 
arithmic scale, where the chain is divided into two equal-sized 
blocks. In contrast to the localized behavior of the participa- 
tion functions in the vacuum sector, the participation func- 
tions in the solitons sector peak at the core of the soliton. 

[IT] we present log Elm of two such sliding blocks, where 
g = 3000, lsoi — 56 and d so i = 281. The positions of the 
solitons are shown for reference. The x coordinate in the 
graph n is the average of the centers of the two blocks A 
and B. Since the double-soliton solution is symmetric, 
the centers of the blocks coincide with the locations of 
the solitons for n — 0. Indeed, the entanglement is max- 
imal at that point. In addition to the main maximum, 
there are secondary maxima at n = ±170, which are ob- 
tained when n coincide with the center of the solitons. 
These maxima correspond to the entanglement between 
the right and left tails of each of the solitons. We note 
that in the vacuum sector two such sliding blocks yield 
practically constant value of long-range entanglement. 

The fact that long-range entanglement has its peak at 
the core of the solitons is also manifested in the struc- 
ture of the Williamson's modes, the so-called participa- 
tion functions [12], [j|: 

zj(n) = Uj(n)vj(n), (23) 

where HaGaUj = XjUj and GaHavj — X^Vj. The 
participation functions qualitatively describe the distri- 
bution of entanglement between two complementary re- 
gions. We split the chain into two equal-sized blocks. In 
the vacuum sector the participation functions are practi- 
cally "localized" [f|. z\ ac , which has the highest contri- 
bution to the entanglement, originates at the boundary 
between the parts. The contribution of the other modes 
(z| ac , z| ac , . . . ) decreases exponentially. The localization 
of Zj ac departs from the boundary as j increases. In 
contrast, in the double-soliton solution, the participation 
functions are not localized and correspond to the entan- 
glement between the solitons. In figure [12] we present 
zf° l , z| o/ , z\ ac and z^ in logarithmic scale, where r\\ 
(the first internal mode) is also plotted for reference. In 
contrast to the participation functions in the vacuum sec- 
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tor, the participation functions in the solitons sector peak 
at the core of the solitons. 



B. Long-range entanglement vs. the separation 
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FIG. 13: (color online). Long range entanglement as a func- 
tion of the solitons separation. The centers of the two blocks 
coincide with the solitons: d — d so i and I — l so i. The 
plots are given for g = 1600, 10 4 . log E LN (g = 1000) of the 
vacuum sector is plotted for reference, where (3(g) = 5.57. 
' ~ yf§ ~ 40, 96 respectively. In the case g — 10 4 , the entan- 
glement of solitons decreases slower than its critical vacuum 
sector correspondence for small values of d/l. In the case 
g — 1600, the entanglement oscillates and is not a decreasing 
monotonic function of the separation for large values of d/l. 

We now look at scaling of long-range entanglement as 
a function of the separation between the solitons. We 
compare long-range soliton-soliton entanglement to the 
corresponding entanglement in the vacuum sector, which 
is equivalent to long-range entanglement in the linear har- 
monic chain. In the critical vacuum sector the scaling of 
the LN for sufficiently large chains is universal and de- 
pends on a single parameter r = d/l. For sufficiently 
large values of r, the logarithmic negativity scales as 
Eln ~ e~ l3cr , where (i c p=> 2.7 [23| is the universal decay 
coefficient of the critical chain. In the noncritical chain 
the logarithmic negativity is not universal. Given two 
blocks with constant length /, their entanglement scales 
as E LN fa e -P(.™,l)d/l for d » I, where (3(m,l) > f3 c and 
m is the mass of the particles. 

In the soliton sector we find unique behavior of long- 
range entanglement, which is qualitatively different from 
that of the vacuum sector. Here we compute entangle- 
ment between two blocks that coincide with the solitons 
(d = d so i) for different values of d so i. The size of the 
blocks I is chosen so that it maximizes von-Neumann en- 
tropy Es (Eq. [IS]) • As expected I does not depend on the 
separation d = d so i and approximately, I ~ y/g. In figure 
[Hwe Plot \ogE LN versus d/l for g = 1600, 10 4 . Both 
values of g describe the noncritical regime. (There is no 



critical regime in which the soliton solutions have inter- 
nal modes.) For reference we also plot \ogE LN (g = 10 4 ) 
for the corresponding vacuum sector with the same value 
of I. Clearly, long-range entanglement in the soliton sec- 
tor is larger than that in the vacuum sector (this is true 
for all values of g). Note also that the scaling of LN in 
the soliton sector differs qualitatively from that in the 
vacuum sector as log Eln is not linearly dependent on 
d/l. For relatively small values of d/l a linear curve ap- 
proximation of log E L ° N (g = 10 4 ) yields (3 — 2.3, which 
is even stronger than the universal scaling of the critical 
vacuum sector 23]. In this regime the solitons have a 
large overlap, corresponding to the negative square mass 
in the region of the solitons. As d/l increases the linear 
curve does not fit and has a somewhat oscillating behav- 
ior. For large values of d/l we see that the decay of the 
entanglement is larger than the corresponding one in the 
vacuum sector. 

For the smaller value of the coupling constant, g = 
1600, we observe even more surprising behavior. The LN 
oscillates around a somewhat constant value for d/lo > 8. 
In this regime long-range entanglement loses one of its 
most fundamental characteristics: it is no longer a mono- 
tonically decreasing function of the distance between the 
solitons! 



C. Toy model for long-range entanglement 
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FIG. 14: (color online). Particle-like entanglement of the soli- 
tons: dependence of log Eln on the ratio between the inter- 
nal modes' eigenfrequencies for g = 10 4 , where 201ogo;i and 
201ogoj2 are also plotted. The three left vertical arrows are 
examples of drastic increase of the gap, which correspondingly 
changes the entanglement. Notice that as the gap is closed 
(d/l > 6), the entanglement decreases significantly. The en- 
tanglement of the toy- model, E(a), is shown for reference. 

We would like to analyze the unique behavior of the 
long-range entanglement in the soliton sector. We discuss 
two factors, each of which dominates in a different regime 
of the coupling constant g. For relatively large coupling, 
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FIG. 15: (color online). The exponential coefficient (3 in the 
asymptotic behavior Eln{<1 » I) ~ e ~' 3d / ! as a function of 
g. In the vacuum sector /3 decreases monotonically to the 
universal value 2.7 (not seen in graph). The lower bound (3 > 
2 is satisfied in the soliton sector for approximately g > 2000. 
For smaller values of g the exponential approximation breaks 
down. 



g > 5000, the particle-like behavior of the solitons is the 
dominant factor. For relatively low coupling, g < 1000, 
it is the topological nature of the first internal mode that 
dictates the entanglement behavior. 

Let us first discuss the particle-like behavior of the 
solitons by considering a toy-model. From figure 2] we 
see that the first two internal modes: 171, 772 are local- 
ized to the solitons and describe qualitatively two groups 
of oscillating particles inside each of the solitons. 771 
and 7/2 are symmetric/ antisymmetric and describes corre- 
lated/anticorrelated oscillations of the two groups respec- 
tively. Therefore, we simplify our model to two coupled 
harmonic oscillators with normal-mode eigenfrequencies 
u)% and u>2, which we identify with the internal modes' 
eigenfrequencies in the solitons sector. Then it can be 
easily shown that the entanglement between the oscilla- 
tors depends on a single symplectic eigenvalue 



A 



where a — LOijio\. Correspondingly, the entanglement 
between the solitons depends mainly on the ratio between 
the two internal modes eigenfrequencies. 

In figure [14] we show the particle-like behavior of the 
two solitons for g = 10 4 . First note that both u>i(d/l) 
and L02{d/l) are not monotonic functions. Each has a cer- 
tain global minimum, which appears due to the finiteness 
of the system. As the gap between the eigenfrequencies 
is closed (d/l > 6), the entanglement decreases signifi- 
cantly. Moreover, the slight deviations from the curve 
appear when there are also deviations in the gap, as in- 
dicated by the vertical arrows. For reference, we present 
the entanglement of our toy mode E(a). Clearly, the 
toy-model manifests much larger entanglement since the 
phonon band is not included, and therefore its decohering 



screening effect is missing. 

This model can be used to obtain a lower bound for 
the decay coefficient (3 in the limit where d >> I. In this 
limit we can use a double-well WKB approximation to 
estimate the eigenfrequencies: 



oc exp 



1/2 



\ -d sol /2 



\p(x)\dx 



(24) 



rpsol 



where p(x) is the classical momentum, which is approx- 
imately 1/y/g. In this case the entanglement can be ap- 
proximated for d >> I as Eg oc e~ 2d /^. Therefore, we 
obtain a lower bound for (3 in the soliton sector, (3 > 2. 
We expect that this bound is valid for large values of 
g, where the particle-like behavior of the solitons is the 
most significant effect, and for large values of d/l so that 
the solitons have a small overlap. However, for large val- 
ues of the coupling constant we cannot study very large 
values of d/l since the system is finite. 

Figure [T5] shows the exponential decay coefficient in 
-pd/l m asymptotic limit d >> I, where the 
exponential approximation becomes accurate. f3(Ef/§) 
is also shown for reference. In the vacuum sector fj de- 
creases monotonically down to the universal value 2.7 
(not seen in graph). For g > 10 4 the double-soliton sec- 
tor has similar behavior. However, for small values of g 
the discussed lower bound is not satisfied, (3 < 2, as we 
enter the regime where u>\ is the most significant factor 
and the exponential approximation breaks down. 



D. The noncritical regime 
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FIG. 16: (color online). Topological nature of wi. logi^Liv 
versus d for g = 731 (deep in the noncritical regime). The 
computed logarithmic negativity is not a monotonically de- 
creasing function of the separation between the solitons. For 
reference we plot 21og(wi) — 20. Note that maximal entan- 
glement is achieved for minimal Wi, which characterizes the 
regime in which g is small. 

Let us now turn to the second factor, which domi- 
nates for relatively low coupling (g < 1000). As described 
in Appendix A the double-soliton solutions depart from 
4>sg (Eq. as we require a finite system and impose 
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FIG. 17: (color online). 1^2/1^1 and logtJi as a function of 
log E LN for three coupling-strengths: g = 29500, 7900, 1000. 



d(f>(±L) j dx 7^ 0. The separation of the solitons d so i is 
given by a parameter fc, where d so i monotonically in- 
creases with k. Stability (wf > 0) regions of the system 
are defined in separate intervals of k. As the values of k 
in an arbitrary interval get closer to the interval bound- 
aries, uj\ becomes closer to zero. Therefore, we obtain an 
oscillatory behavior of lo\ with the separation d so i , which 
strengthens as discreteness increases (g decreases). In 
figure [TBI we present log Eln as a function of the sep- 
aration between the solitons d so i, where g = 731. For 
reference we present on the same graph 2 log(ci; 1 ) — 20 so 
that both plots are more or less on the same magnitude. 
Note that maximal entanglement is achieved for mini- 
mal wi . Clearly, in this regime long-range entanglement 
is no longer a monotonically decreasing function of the 
distance between the solitons. 

In figure [17] we show both W2/^i and logwi as a 
function of log E^n for three coupling values g = 
29500, 7900, 1000. In the close to critical regime (g = 
29500) the particle-like behavior of the solitons is a dom- 
inant factor as it is a monotonic function of the entan- 
glement, where the first internal mode has no significant 
effect (its minimal value does not correspond to maxi- 
mal entanglement). In the deep noncritical regime the 
particle like behavior has no significant effect as a->2 /^i is 
approximately constant (slightly above 1), and the dom- 
inant factor is w%. In the case g = 7900, we see that 
in general, the entanglement increases with wilui\ and 
decreases with uj\. However, there is an intermediate re- 
gion, where the entanglement is not a monotonic function 
of both factors. 



V. THE WEAK-COUPLING LIMIT 

We turn now to the weak-coupling limit of the FK 
model, which is characterized by small values of g. Here 
the entanglement has far less dramatic effects. No long- 
range entanglement has been detected between two kinks. 
In addition, the von Neumann entropy of an increasing 
block inside the kink loses its localization as the cou- 
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FIG. 18: (color online). Es as a function of / for small values 
of g. 



pling constant becomes too low, as seen in figure [TBI For 
0.25 > g < 2.25 there are two internal modes in the soli- 
ton sector of the FK model [lj| , where the second internal 
mode has an eigenfrequency close to one. However, the 
additional internal mode has no significant contribution 
to the entanglement of the kinks. 



VI. SOLITONS AS CARRIERS OF QUANTUM 
INFORMATION 



In this section we suggest and investigate the possible 
use of solitons as carriers of quantum information. We 
observed in the previous sections that the entanglement 
properties of solitons in their ground state are strongly 
related to the presence of internal localized modes. These 
are (usually) the lowest discrete modes in the spectrum of 
the normal modes. The internal modes qualitatively de- 
scribe collective or vibrational perturbations of the clas- 
sical soliton. In our approximation method the soliton is 
handled as a non-perturbative classical object while the 
internal modes are linear perturbations around it. The 
quantized internal modes then correspond to "phonon" 
bound states. 

Our suggestion therefore is to use the localized inter- 
nal modes' degrees of freedom as carriers of quantum 
information. Since the modes are "attached" to the soli- 
ton, by moving the soliton under suitable adiabatic con- 
ditions, information may be transported without being 
disturbed or mixed with the phonon band. It is also pos- 
sible to formulate multi-soliton processes which create 
entanglement, and therefore manipulate the information 
carried by each of the solitons. 

To our knowledge this idea has not been investigated. 
Therefore we provide a preliminary study of this possibil- 
ity without elaborating on a particular implementation. 
In our model information will be coded by exciting the in- 
ternal mode(s). A fundamental requirement in quantum 
information processing is the possibility of local opera- 
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tions. However, the internal modes are collective degrees 
of freedom of all the particles in the system. To resolve 
this difficulty one has to show that a local addressing of 
the particles in the core of the soliton is sufficient in or- 
der to manipulate the internal mode(s) with high fidelity. 
Such a process seems possible since the internal modes 
are localized. 

To support the above requirement we present two re- 
sults regarding the single soliton and double-soliton re- 
spectively: 

(1) Suppose we perform two-mode-squeezing of the in- 
ternal mode with some arbitrary external degrees of free- 
dom, Q. Then the entanglement inserted by squeezing is 
saturated by the entanglement of local degrees of freedom 
inside the soliton. 

(2) Suppose we perform two-mode squeezing of two 
new collective modes that are linear combinations of the 
internal modes, such that the new modes are confined to 
the left and right sides of the chain respectively. Then 
the entanglement inserted by squeezing is saturated by 
local degrees of freedom inside each of the two solitons. 

Both of these results are crucial for any quantum infor- 
mation applications, as otherwise operations cannot be 
implemented locally. 



A. Single soliton 



where Eq. [10] can be rewritten as <j> = rjcj) and it = r/7r 
in matrix form. The degrees of freedom of the internal 
mode are ffi). The ground state of the system is 
clearly 



1-0. 



1=1.2,. ..,N 



i=i 



(26) 



Two- mode squeezing S(r) of the internal mode and Q is 
defined by 



S(r) 



^r(a{a { Q -a iaQ ) 



(27) 



where r is the squeezing parameter. The squeezing acts 
on the internal mode and Q directly through: 



4>i -> (e +, '0i + e- r x Q )/V2, 
7r 1 -»(e-'-7ri + e +r po)/V2, 
X Q ( e+r 4>i - e~ r x Q )/V2, 
PQ ->■ (e 7Ti - e T pq)/V2, 

We can plug S(r) into the covariance matrix 



(28) 



M -> Rc (r}2 1 SMS -X J7 2 ) , 



(29) 



where 772 = ?7©1©?7©1 and 




■E S (A)-E S (B) 

E s(°) 



400 600 
I 



800 1000 



FIG. 19: (color online). Entangling a soliton and an external 
oscillator Q by two-mode squeezing with r — 0.99. Shown 
E S {A) - E S {B) < E D (A) as a function of I for g = 1000. 
Es{Q) is shown for reference. Note that the distillable entan- 
glement saturates the inserted entanglement for quite small 
values of I. 

Let us start with the first requirement. We write the 
normal modes of the system as: 



(25) 



Til = -l< 



1 /If 11 

M = -Diag — , — , . . . , — , - ,Wi,cj 2 , ■ ■ ■ ,wn,wq 

2 Vwi w 2 wn 0JQ 



is the diagonal vacuum covariance matrix in the normal 
modes. 
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FIG. 20: (color online). Saturation of the hashing inequal- 
ity - schematic presentation. At t = we perform two-mode 
squeezing of x an d Q so that x decay significantly outside 
the squeezing block W and we measure the entanglement be- 
tween A and Q such that A » W . In a simplified mode- 
wise decomposition of A and BuQ we treat only two modes: 
| 1a, Is), in which the participation function is localized at the 
boundary and \2a, 2b), in which the participation function is 
localized inside W. 
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We show that for sufficiently large block A, which in- 
cludes the soliton, the distillable entanglement E D be- 
tween local degrees of freedom inside A and Q saturates 
the inserted entanglement via squeezing: 

Es = cosh 2 r ln(cosh 2 r) — sinh 2 r ln(sinh 2 r). (30) 

From hashing inequality [24| . 

S D (A,Q)>max{0,E s (A)-E s (AuB),E s (A)-E s (AuQ)}, 

(31) 

where B complements A in the chain so that 

E S (AUB) = E S (Q), 

E S (AUQ) = E S {B). 1 ' 

In figure [TO1 we plot Eg(A) — Es(B) as a function of I 
for g — 1000. Eg{Q) is shown for reference. Note that 
indeed E D (A,Q) saturates Es(Q) for I > 450. 

We would like to provide a qualitative explanation of 
the saturation obtained in figure [T91 Note that the dis- 
cussed saturation does not depend on the existence of an 
internal mode, but on the shape of the squeezing. The 
internal modes are required for the entanglement to stay 
in the squeezing window, but at t = they do not have 
a significant role, where an arbitrary collective mode x, 
which has the shape of a bound state, will be saturated 
by local modes. We therefore turn to analyze the vacuum 
sector. We specify two blocks in the chain as can be seen 
in figure[2ni the squeezing block W and the entanglement 
measuring block A, where A >> W. We measure the en- 
tanglement between Q and A. In addition, let us assume 
that x decays significantly outside W. Saturation of the 
hashing inequality means: Es(A) — Es{B) — > Es{Q). 

As the whole state is Gaussian we use the mode-wise 
decomposition That is, the state is described by a 

product of composed collective modes from A and BUQ: 

IVO = Y[\i>i)A,(BUQ)- (33) 

Now before squeezing, the most significant contributions 
to Es(A)(= Es(B)) come from collective modes that are 
localized at the boundaries between A and B. The lo- 
calization of less significant modes tends to the center of 
A. Therefore, modes that are well concentrated inside 
W contribute very little to the initial entanglement. The 
hashing inequality is saturated if after squeezing, Es{B) 
does not change much, that is, almost all entanglement 
that is inserted by squeezing is added to Es(A). 

To see that this is the case consider the following sim- 
plified model. Let us assume that before squeezing the 
state can be expressed as a product of two composite 
states: 

\f) = |U,1b)®|2 A) 2 b ) (34) 

where |Ea,1b) is a strongly entangled mode, such that 
| 1a) is localized at the boundaries of A (outside W) and 
|2a,2b) is a weakly entangled mode, such that \2a) is 



localized inside W . Squeezing entangles \2a) to Q but 
it only weakly entangles \1a) to Q. The entanglement 
between B and Q is therefore modified by two negligible 
contributions: one from weak squeezing (Is) and the 
second from weak entanglement in the first place (2b). 

B. Double-soliton 
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FIG. 21: (color online). Double soliton solution. r/± — 
-^(Vi ± m) for 9 = 500 and g = 10 4 , where u>i(c? = 500) = 
0.0011, ui 2 (g = 500) = 0.0012, ui(g = 10000) = 0.0832 and 
tJ2(g ~ 10000) = 0.3075. Note that for almost degenerate 
modes the transformed collective modes r/± are well sepa- 
rated, while for nondegenerate modes the collective modes 
overlap and the squeezing is not optimal. 
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FIG. 22: (color online). Double soliton solution. Two- 
mode squeezing of r]± with r = 2. The graph shows 
Es(Ai) - Es{Ai U A 2 ) < E D (Ai,A 2 ) as a function of I for 
g = 500, k = 0.9415 and d(l) = 331 — I (the separation be- 
tween the blocks decreases as their size increases). E(±) is 
plotted for reference. Note that the distillable entanglement 
actually exceeds the inserted entanglement by squeezing, as 
for large values of I, d(l) decreases to very small numbers, for 
which the vacuum entanglement has a non-negligible contri- 
bution. 
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We turn now to the second result regarding double- 
soliton solutions. Now the two internal modes are sym- 
metric and antisymmetric respectively. They are not con- 
fined to a single soliton each. We want to create entan- 
glement between the solitons by squeezing. Therefore, 
we should use collective modes, where each is confined to 
a single soliton. Let us define the linear transformation: 

v± = -^(vi±m)- (35) 

In the case the internal modes are almost degenerate, 
rj± are localized well to each of the solitons. Note that 
this condition corresponds to poor long-range vacuum 
entanglement of the solitons, as described in Section 4. 
Practically degenerate internal modes appear only deep 
in the noncritical regime of the doublc-soliton sector (g < 
1000), as shown in figure |2"TI where we plot r/± for g — 500 
and g = 10 4 . 

We therefore assume that g = 500. In correspondence 
with the single soliton squeezing, we define 

S(r) = e »-(«V- -«+»-) (36) 

as optimal squeezing of the two solitons. We would like to 
show that the inserted entanglement by squeezing E(±) 
(Eq. |3"0|) is saturated by local degrees of freedom inside 
the two solitons: A\ and A 2 with the same size. From 
symmetry, Es{A\) — Es{A 2 ). Hashing inequality im- 
plies that, 

E D (A U A 2 )>E S (A 1 )-E S (A 1 UA 2 ). (37) 

In figure [221 we plot Es(Ai) — Es(A± U A 2 ) as a function 
of blocks size I, where the center of the blocks coincide 
with the center of the solitons. E(±) is plotted for ref- 
erence. Note that the distillable entanglement saturates 
the inserted entanglement for I > 250. Therefore, lo- 
cal degrees of freedom saturate E(±). Interestingly, the 
distillable entanglement actually exceeds the inserted en- 
tanglement through squeezing. This is explained by the 
fact that for large blocks the distance between them be- 
comes very small, such that the vacuum entanglement 
has a non-negligible contribution. 

C. Implementations 

Finally, based on the above analysis, we would like to 
briefly discuss two possible applications: entanglement 
transportation and entanglement manipulation through 
a tunneling gate. In the entanglement transportation 
task we assume there is an additional external degree of 
freedom, Q. Initially the soliton is located in one side 
of the system, next to, say, Alice. Alice entangles Q to 
some local degrees of freedom inside the soliton. Then she 
classically transport the soliton to the other side, next 
to Bob. Finally, Bob possesses local degrees of freedom 
inside the soliton with which he can, e.g., perform ex- 
periments that prove nonlocal correlations with Q. In 



order to keep the state Gaussian, the entanglement is re- 
alized by two-mode squeezing. In addition, we assume 
slow transportation so that the quantum modes are per- 
turbed only adiabatically. 

In the second implementation we use the double- 
soliton sector. Quantum information is tunneled between 
the solitons by bringing them close to each other for a fi- 
nite time interval and then separate them apart. This 
application is much more complicated than the first one, 
as we require that the phonon band is not excited during 
the tunneling phase, while the internal modes do. 



VII. CONCLUSIONS 

We conclude with the main results in this work and 
proposals for further study. We have shown that the 
entanglement of solitons manifests particle-like behav- 
ior as they are characterized by localization of entangle- 
ment and long-range entanglement. The von-Neumann 
entropy of solitons mixes both critical and noncritical 
behaviors, where it increases logarithmically inside the 
core of the soliton, reaches a maximal value and then 
decreases and saturates the constant value that corre- 
sponds to the massive vacuum sector. Interestingly, the 
increase of the entanglement inside the core is faster than 
that in a universal critical field. In addition, we have 
shown that two solitons manifest long-range entangle- 
ment, which may not decrease with their separation in 
the noncritical regime. Near the critical regime the en- 
tanglement decreases slower than the corresponding uni- 
versal critical field. 

Our quantum model is based on linear perturbations 
around a nonlinear classic solution. Higher order cor- 
rections [25j are important for the treatment of the zero 
mode in the infinite system. In this case one obtains 
non-linear terms which couple the solitons' collective co- 
ordinate with the phonon modes. It would be interesting 
to explore classical configurations of two non-static soli- 
tons which can produce entanglement through their mu- 
tual phonon modes. Moreover, certain models, such as 
the (j) 4 , manifest in addition to the zero mode a nonzero 
internal mode which is associated with local shape de- 
formations. It would be interesting to study mechanisms 
for entanglement of the internal modes in soliton-soliton 
scattering [26| and in soliton- impurity scattering (27| . 

We suggest that solitons' internal modes be used as 
carriers of quantum information. Since the modes are 
"attached" to the soliton, by moving the soliton un- 
der suitable adiabatic conditions, information may be 
transported without being disturbed or mixed with the 
phonon band. Our proposal differs from previously 
suggested methods, which employ a collective displace- 
ment of the soliton using long Josephson junctions [28| . 
Realizations of the present suggestion can possibly be 
achieved either in discrete systems such as trapped ions 
or in an effectively one-dimensional Bose-Einstein con- 
densate, which manifests soliton solutions [29|, [3(j. Re- 
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cently, a possible realization of an FK-like model in the 
ion trap was discussed in [31( . The utilization of inter- 
nal modes in trapped ion systems will be discussed in a 
future work 1331. 
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APPENDIX A: SOLUTIONS OF THE 
SINE-GORDON EQUATION FOR FINITE 
SYSTEMS 

For a comprehensive analysis see [33]. The finite 
sine-Gordon system describes the magnetic flux in long 
Josephson junctions with constant, homogeneous exter- 
nal magnetic field H perpendicular to the barrier and 
bias current line density /. We use dimcnsionlcss units: 

x^x/Xj, H^H/H S , I -»///., 

where Aj is Josephson length and 

H S = I S = [e(2A L +f)A,/] _1 

is the superheating field of the vortex- free Meissner state, 
where t is the thickness of the barrier and Al is London 
penetration depth. 

In the finite system the boundary condition we impose 
on the magnetic flux <j> is: 



J(±L) = ±I- 
ax 



2H 



where 2L is the length of the barrier in Josephson's length 
units. Imposing stationary condition on sine-Gordon 
equation, Eq. [4] reduces to 

d 2 (f>/dx 2 = sin</>, 

where for zero transport currents the solutions are: 

<f>Jx) = tt(ct - 1) + 2am(^ + K(k),k), a = 0, 2, 4, . . . , 

k 

x 

4> (x) — tt<7 + 2am( — , k), a = 1,3,5, ... , 
k 

(A.2) 

where am is the Jacobi elliptic amplitude and K is the 
complete elliptic integral of the first kind. The subscripts 
e (even) and o (odd) refer to the number of soliton (vor- 
tices) in the solutions. 

These solutions are stable only for specific conditions 
on k, L and H . Stability regions in terms of k and H 
for Eq. IA.2I are established from the following: first, the 
roots of the equation 

aK{k)=L, a = 1,2,..., 



form an infinite decreasing sequence of bifurcation points, 
where a is the analog of the topological index a in the 
infinite single soliton sector. 

k = k a e 1= (0,1] = u~ J ff , 



Itx — (ka+i i k a ] , 

where a = 1, 2, . . . The stability regions for the solutions 
(j) = 4> e are given by the intervals Ii m , whereas the sta- 
bility regions for the solutions cj> = <j) are given by the 
intervals l2 m +i, where m = 0, 1, 2, . . . Dependence of k 
on H is given by 



i-fc 2 

kH 



-, a = 2m, 



(L 



(A.3) 



dn\^—,kj=kH, a = 2m + l, 



where 



—dn(u, k) = —k ■ cn(u, k)sn(u, k) 
ou 

where cn(u, k) and sn(u, k) are correspondingly the el- 
liptic cosine and sine functions. The stability regions in 
terms of the field H take the form: 



< H < H , a = 



Hl_ x - 1 < H < H a , a =1,2,... 



(A.4) 



(A.l) with H a implicitly determined by 



(a+l)K^—j=H a L. 

Turning to the quantization of the perturbation, one 
has to solve Eq. [H] 

— ?/' + [2ksn 2 (x + Xq, k) — l] r\ = u> 2 r] (A. 5) 

where xo is the phase of the center of the junction. Eq. 
IA.5l is also known as Lame equation. The boundary con- 
ditions for rj are rf\ x —±L = |34j ]. 

Finally we would like to note that if bias current were 
introduced, the locations of the solitons in the presented 
solutions could be modified to a nonsymmetric config- 
uration. In addition, note that g does not appear in 
the equations. Recall that in the infinite single soliton 
solution the width of the soliton is ^fg. In the contin- 
uous finite-system case one can not change arbitrarily 
the width of the vortices, as they scale proportionally 
with the system size. Throughout the analysis we choose 
L = 8 as the sampled analytical solution we start with 
in the minimizing algorithm. 
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